function F = sys(x,eC,ell,theta,sigma,pi,eta,i,M,A,B)
% endogenous asset specialization
% partial equilibrium, given eC: 0 < eC < 1

q_star = 1;

q0A = x(1);
q1A = x(2);
q0B = x(3);
q1B = x(4);
qBA = x(5);
eN  = x(6);
dBA = x(7);

a_CA_n = eN*(1-ell)/(eC*ell+eN*(1-ell))^(1-eta);
a_CB_n = (1-eN)*(1-ell)/((1-eC)*ell+(1-eN)*(1-ell))^(1-eta);
a_NA_n = eC*ell/(eC*ell+eN*(1-ell))^(1-eta);
a_NB_n = (1-eC)*ell/((1-eC)*ell+(1-eN)*(1-ell))^(1-eta);
a_CA_d = 1-ell;
a_CB_d = 0;
a_NA_d = ell;
a_NB_d = 0;

F(1) = (qBA - q0B) * (ell*(pi*a_CA_n+(1-pi)*a_CA_d)*theta/w(q1A,sigma,theta)*(mu(q1A,sigma)-1) ...
                     -ell*(1-pi)*a_CA_d*theta/w(qBA,sigma,theta)*(mu(qBA,sigma)-1)); 

F(2) = - i + ell*(1-(pi*a_CA_n+(1-pi)*a_CA_d)*theta/w(q1A,sigma,theta))*(mu(q0A,sigma)-1) ...
           + ell*(pi*a_CA_n+(1-pi)*a_CA_d)*theta/w(q1A,sigma,theta)*(mu(q1A,sigma)-1);
F(3) = - i + ell*(1-pi*a_CB_n*theta/w(q1B,sigma,theta)-(1-pi)*a_CA_d*theta/w(qBA,sigma,theta))*(mu(q0B,sigma)-1) ...
           + ell*pi*a_CB_n*theta/w(q1B,sigma,theta)*(mu(q1B,sigma)-1) + ell*(1-pi)*a_CA_d*theta/w(qBA,sigma,theta)*(mu(qBA,sigma)-1);

F(4) = - q1A + min(q_star, q0A + ((eC*q0A+(1-eC)*q0B)/M*(A-(1-eC)*dBA)/eC - (1-theta)*(u(q1A,sigma)-u(q0A,sigma)))/theta);
F(5) = - q1B + min(q_star, q0B + ((eC*q0A+(1-eC)*q0B)/M*B/(1-eC)          - (1-theta)*(u(q1B,sigma)-u(q0B,sigma)))/theta);
F(6) = - qBA + min(q_star, q0B + ((eC*q0A+(1-eC)*q0B)/M*dBA               - (1-theta)*(u(qBA,sigma)-u(q0B,sigma)))/theta);

S_NA   = (1-theta)*(u(q1A,sigma)-u(q0A,sigma)-q1A+q0A);
S_NB_n = (1-theta)*(u(q1B,sigma)-u(q0B,sigma)-q1B+q0B);

F(7) = a_NA_n*S_NA - a_NB_n*S_NB_n;

%=========================================================================
% subfunctions
%-------------------------------------------------------------------------
function u = u(q,sigma)
% CRRA utility: u(q)
if sigma < 1
    u = q.^(1-sigma)./(1-sigma);
elseif sigma == 1
    u = log(q);
end

function mu = mu(q,sigma)
% marginal utility: u'(q)
mu = q.^(-sigma);

function w = w(q,sigma,theta)
% w function related to effective bargaining power
w = theta + (1-theta)*mu(q,sigma);
